Atomistic mechanism of physical ageing in glassy materials 
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Abstract. - Using molecular simulations, we identify microscopic relaxation events of individual 
particles in ageing structural glasses, and determine the full distribution of relaxation times. We 
find that the memory of the waiting time elapsed since the quench extends only up to the 
first relaxation event, while the distribution of all subsequent relaxation times (persistence times) 
follows a power law completely independent of history. Our results are in remarkable agreement 
with the well known phenomenological trap model of ageing. A continuous time random walk 
(CTRW) parametrized with the atomistic distributions captures the entire bulk diffusion behavior 
and explains the apparent scaling of the relaxation dynamics with tw during ageing, as well as 
observed deviations from perfect scaling. 



The dynamics of glass forming materials is characterized 
by slow, spatially heterogeneous and temporally intermit- 
tent processes [1-6] . Long periods of caged particle vibra- 
tions are interrupted by rapid, localized structural relax- 
ations. Dynamical heterogeneity near the glass transition 
leads to stretched exponential relaxation of correlation and 
response functions, as well as subdiffusive, non-Gaussian 
transport. In the glassy state, configurational degrees of 
freedom continue to evolve slowly. This phenomenon is 
called physical ageing [7] and is a defining and universal 
. feature of nonequilibrium glassy dynamics [8] . Structural 
■ relaxations become increasingly sluggish with the waiting 
. time tw elapsed since vitrification, and almost all material 
properties change as a result. 

Ageing often manifests itself as a simple rescaling of the 
slow, configurational part of the global correlation func- 
tions C{t,tni) = Cyib{t) + Cconf{t/t'i,), where ^ is the 
ageing exponent. This scaling behaviour has been ob- 
served in a wide range of disordered materials [1] and is 
particularly pronounced in polymer glasses [7], colloidal 
suspensions [9,10], and physical gels [11]. Recent simu- 
lations of structural glasses show that scaling also applies 
approximately to the distribution of local correlation func- 
tions [12] and the distribution of particle displacements 
(van Hove function) [13], which superimpose when plot- 
ted at times of equal mean. Scaling of local correlations 
can be proved to hold exactly for certain spin glass mod- 
els whose dynamical action obeys time reparametrization 



invariance [14]. Despite the apparently universal nature 
of this phenomenon [15], the physical mechanism behind 
it remains elusive. 

In this Letter, we directly observe the atomistic dynam- 
ics underlying physical ageing on a single particle level 
through molecular dynamics simulations of two model sys- 
tems for atomic and polymeric glasses. We consider both 
the well-known amorphous 80/20 binary Lennard- Jones 
mixture [16], as well as a bead-spring model for polymer 
glasses [17,18], where particles interact via the Lennard- 
Jones potential and are coupled together by stiff springs 
to form chains of length 10. Twenty thousand particles 
each are simulated in a periodic simulation box. For the 
polymer (80/20) model, the initial ensemble is equilibrated 
at T = 1.2 (4.0), and then quenched rapidly to T = 0.25 
(0.33) to form a glass. The system is aged at constant pres- 
sure for various waiting times tw and then the dynamics 
are followed as a function of the measurement time t. All 
results are quoted in Lennard- Jones units. 

In both models, the effect of ageing can be immediately 
observed through the one-dimensional mean squared dis- 
placement {Ax'^{t,t^)) shown in Fig. 1, which exhibits 
three regimes: a short time vibrational regime precedes 
a long plateau characterized by "caged" motion, followed 
by a cage escape or a-relaxation regime. The cage escape 
time increases with increasing tw, but { Ax {t,tw)'^) -curves 
for different waiting times can be superimposed in the cage 
escape regime by rescaling time with t!^. At long times 
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Fig. 1; Mean squared displacement for (a) the binary Lennard- 
Jones mixture and (b) the polymer model for waiting times 
= 750, 2250, 7500, 22500, and 75000 (left to right). From 
superposition in the cage escape regime we deduce fi = 0.51 
for the Lennard- Jones mixture and /i — 0.62 for the polymer 
model. In the polymer model the terminal slope is smaller due 
to Rouse dynamics. 

t 3> tuj, scaling of the mean squared displacement fails 
for both models, and curves for the different ages merge. 
Experimental creep compliance curves show the same be- 
havior, which is a well-known consequence of subageing 

i^,<^) [7]. 

We proceed to analyze the intermittent dynamics in 
these systems by decomposing single particle trajectories 
into discrete hopping events. To investigate the hopping 
dynamics, a subset of five thousand particle trajectories 
are recorded. A running average of the atom position and 
its standard deviation is found for each atom over a small 
time window of 400 (40 individual snapshots of the par- 
ticle positions) [19]. A hop is identified when the stan- 
dard deviation a is greater than a threshold of 2((t) [20]. 
This method differs from the usual technique of identify- 
ing hops via a threshold in the hop distance [5,21,22], and 
effectively captures hops at all length scales through their 
increased activity during the relaxation event. While the 
correlation factor between hops is nearly zero, approxi- 
mately 10% of hops are either forward or backward corre- 
lated. Removing these hops does not affect the results. 

Fig. 2 shows a typical atomic trajectory with the hops 
identified using this criterion. The initial hop times ti, 
the times between all subsequent hops r, also called the 
persistence time, and the displacements dx were binned 
to form distribution functions p(ii), p{t), and p{dx). The 
initial hop time distribution is shown in Fig. 3(a) for the 
Lennard- Jones mixture and in (d) for the polymer. For 
all waiting times, it appears to take the form of two power 
laws. The exponents are waiting time independent, but 
the crossover time between the two regimes shifts towards 
longer times with increasing i^,. 

Figures 3(b) and (e) show that the distributions p{t) 
also decay as a power law: p(r) ~ t--1-20 (Lennard- Jones 
mixture) and p(t) ^ r^^-^^ (polymer). By contrast, p{t) 
is independent of and shows no ageing effects. The one- 
dimensional hop displacement distribution (Fig. 3(c) and 
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Fig. 2; Displacement (red/gray) and standard deviation 
(black) from the running average of a typical atomic trajec- 
tory. The threshold in the standard deviation used for identi- 
fying hops is indicated by a horizontal dashed line. The times 
where hops are identified are indicated by vertical dashed lines. 

(f)) is sharply peaked at zero displacement, and is also in- 
dependent of waiting time. This is perhaps not surprising, 
as it has been shown that the density and the dynamical 
correlation lengths vary little during ageing [6]. p{dx) has 
a purely exponential form for the mixture, whereas in the 
polymer model it decays more rapidly in the tails due to 
chain connectivity effects. Both the first hop time distri- 
bution, and the displacement distribution are in qualita- 
tive agreement with results obtained by analyzing hops 
between adjacent inherent structures [23-25]. 

The measured distributions have much in common with 
the well-known phenomenological trap model of ageing 
[26]. In the trap model, as in our results, p{t) is con- 
stant and p(ti ) continuously evolves towards longer times 
in the ageing regime. In fact, even the form of our mea- 
sured distributions is similar. The distributions p{t) and 
p{ti) decay as weak power laws at long times with an 
exponent —(1 -t- x), where a; < 1 is the noise tempera- 
ture which governs activated hopping between states [27] . 
In this formulation, x is either a constant (exponential 
density of states), or more realistically [23-25], a slowly 
varying function of time (Gaussian density of states) [26] . 
We do not observe any curvature in p{t), but cannot rule 
out its existence at longer timescales. While our results 
clearly show subageing behavior, the ageing exponent in 
the trap model is strictly unity [26]. However, this holds 
true only in the limit ^ oo. We have observed that 
depending on the quench conditions, there can be a long 
lived transient regime where correlation functions approx- 
imately scale with an ageing exponent fi < 1. 

To confirm that our analysis captures all important 
physics, we recompute the ageing dynamics through 
stochastic realizations of a continuous time random walk 
(CTRW) using the measured hop distributions p(<i,t^), 
p(t) and p{dx) as input. For the polymer model, a pure 
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Fig. 3; Distributions of first hop time ti, persistence time r, and displacement dx for the binary Lennard- Jones mixture (a-c) 
and the polymer model (d-f) for waiting times t„ = 750 (Q), 2250 (□), 7500 (A), 22500 (0), and 75000 (<). Dotted lines 
indicate power laws and have slopes of -0.55 and -1.1 (a), -1.20 (b), -0.65 and -1.1 (d), and -1.23 (e). 



random walk is insufficient as the chain connectivity in- 
duces subdiffusive Rouse dynamics. The mean squared 
displacement of the atoms as a function of the number 
of hops n is found to be (Ax^) n^^^ ((Ax^) ^ n for 
the binary mixture). In the CTRW, the Rouse behaviour 
is accounted for by choosing a new particle position from 
a Gaussian distribution widening as n^^^ at every step. 
This procedure creates the proper displacement statistics. 
In both systems, the initial displacement is chosen such 
that the mean squared displacement coincides with the 
MD data at the shortest measured relaxation time. This 
is slightly higher than the actual vibrational amplitude, 
but accounts for relaxations which occur at times down 
to the vibrational timescale, but could not be captured 
otherwise. 

The full distribution of displacements (van Hove func- 
tion) 

Gs{x,t,t^) = {8(x - \xi(t,t^) - Xi(0,tu,)|)) (1) 

as well the as the mean squared displacements are com- 
pared with MD data in Fig. 4 for the polymer model. As 
time progresses, a Gaussian caged peak at small displace- 
ments loses particles to widening exponential tails. The 
origin of the exponential tails has recently generated in- 
tense interest, as they appear to be a universal feature 
of structural glasses [22,28-31]. In our calculation, they 
arise naturally from the wide distribution of persistence 



times. The CTRW reproduces the true ageing dynamics 
exceedingly well, with no fit parameters. The same con- 
clusion holds for the binary mixture. As in the trap model, 
the power law form of p(t) is sufficient to generate age- 
ing. Distributions p{ti,t^) evolved via the CTRW from 
the shortest waiting time alone closely track measured dis- 
tributions. Ageing emerges here as a self-generating, dy- 
namical phenomenon. 

Our results clearly show that there are two relevant 
timescales in the relaxation dynamics of the glass: the 
first hop time which experiences ageing, and the subse- 
quent hop times which are tuj independent. This would 
seem to refute the hypothesis that the dynamics are uni- 
versally reparametrized with waiting time as has been sug- 
gested for spin glasses [12], and may serve to explain some 
deviations from local dynamical scaling that have been ob- 
served. In ref. [12], the authors note that the distribution 
of local correlation functions for a binary Lennard- Jones 
glass superimposes to first order at times where the global 
correlation function is equal. However, good collapse was 
found only for large and small global correlation with sig- 
nificant deviations for C = 0.3 and C = 0.5. 

We further investigate the local dynamical scaling by 
using the CTRW model to calculate the van Hove distri- 
bution for three waiting times and three different values of 
the mean in Fig. 5. Note that here p(ti,tu,) is generated 
by evolving the system directly in the ageing CTRW. The 
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Fig. 4: Van Hove function for the polymer model at — 
75000, and times of 400, 4000, 40000, 400000 (left to right). 
Solid lines are molecular dynamics data and symbols are the 
results of the CTRW parameterized with the measured distri- 
butions. Inset shows the agreement between the mean squared 
displacement from MD (lines) and CTRW (symbols) at ti^ — 
750, 7500, and 75000 (left to right). 



Fig. 5: Van Hove functions from ageing CTRWs at three 
different waiting times 750 (solid), 7500 (dotted), and 75000 
(dashed) plotted at times t where (Aa;^) = 0.03, 0.3, and 3 as 
indicated by horizontal lines in the inset. For these values, the 
correlation function C{t,tw) = {exp{7.2i[x{tw +t) — x{tjn)])) ~ 
0.65, 0.3, and 0.1. Data is generated from the distributions 



p{ti,tw = 0) = 



and ' 



superposition of the van Hove functions displays the same 
behavior reported in ref. [12]. At short times (large C), the 
van Hove functions closely superimpose at times of equal 
mean. In this regime, few atoms have hopped and their 
mean squared displacement and van Hove function are 
dominated by the distribution p(ti, t^), which scales with 
tw For the mid-range of the correlation functions, the 
overlap of the van Hove functions becomes weaker. The 
long waiting time distributions have a sharper trapped 
peak and wider tails than the short function. This is 
due to the fact that the width of the tails is controlled by 
p{t) , whereas the height of the caged peak is controlled by 
p{ii,tw)- At very long times t^t^ (small C), the mean 
number of hops per atom is much greater than one, and 
the waiting time independent behavior begins to dominate 
the dynamics. In this case, the mean squared displacement 
curves merge and there is trivial superposition of the van 
Hove functions. 

To summarize, we have presented a complete charac- 
terization of particle diffusion in ageing structural glasses. 
We find that on a microscopic level, time is not simply 
reparameterized during ageing. Only the first hop time de- 
pends on the waiting time, and the particle subsequently 
"forgets" its age. The evolution of the first hop time dis- 
tributions is remarkably close to that assumed in the trap 
model with annealed disorder. Although this model is a 
rather abstract picture of glassy dynamics, our results in- 
dicate that it may be closer to physical reality on a single 
particle level than previously thought. An explicit map- 
ping of the MD data onto a CTRW with measured dis- 
tributions of hop times and displacements is able to com- 
pletely reproduce all dynamical features of the relaxing 
glass, and provides an explanation for the scaling with 



and small deviations thereof. 
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